#include "exportVTK.h"

void exportSolutionToVTK(TMesh* pTMesh_data, string filename){
	// open file
	ofstream fid;
	fid.open(filename.c_str());

	// check if it's OK, otherwise terminate program
	if ( !fid.is_open() ){
		char msg[512]; sprintf(msg,"File '%s' could not be opened or it does not exist.\n\tCheck if directory was typed correctly.\n",filename.c_str());
		throw Exception(__LINE__,__FILE__,msg);
	}

	int m,n;
	pTMesh_data->allVertices.getSize(m,n);
	// print data to file
	fid << "# vtk DataFile Version 2.0\n";
	fid << "Two phases flow simulation\n";
	fid << "ASCII\nDATASET UNSTRUCTURED_GRID\n";
	fid << "POINTS " << m << " float\n";

	printVerticesCoordenates(fid,pTMesh_data);
	printElementConnectivities(fid,pTMesh_data);
	printCellTypeList(fid,pTMesh_data);
	printPressure(fid,pTMesh_data);
	printSaturation(fid,pTMesh_data);
//	printVelocity(fid,pTMesh_data);
//	printPressureGradient(fid,pTMesh_data);
//	printSaturationGradient(fid,pTMesh_data);
	fid.close();
}

// print vertex coordenates and transfer to each one all computed values to be printed later
void printVerticesCoordenates(ofstream &fid, TMesh* pTMesh_data){
	int nrows, ncols;
	pTMesh_data->allVertices.getSize(nrows,ncols);
	for (int j=0; j<nrows; j++){
		VertexData* vdata = pTMesh_data->allVertices.getValue(j,0);
		fid << vdata->x << " " << vdata->y << " " << vdata->z << endl;
	}
}

// print elements connectivities
void printElementConnectivities(ofstream &fid, TMesh* pTMesh_data){
	const int dim = 1;;
	const int numElements = pTMesh_data->getNumLocalEdges();
	fid << "\nCELLS " << numElements << " " << (dim+2)*numElements << endl;
	int m,n,i;
	pTMesh_data->allEdges.getSize(m,n);
	for(i=0;i<m;i++){
		VertexData* v_0 = pTMesh_data->allEdges.getValue(i,0);
		VertexData* v_1 = pTMesh_data->allEdges.getValue(i,1);
		fid << 2 << " " << v_0->ID_local << " " << v_1->ID_local << endl;
	}
}

void printCellTypeList(ofstream &fid, TMesh* pTMesh_data){
	const int numElements = pTMesh_data->getNumLocalEdges();
	fid << "\nCELL_TYPES " << numElements << endl;
	for(int i=0; i<numElements; i++){
		fid << 3 << endl;
	}
}

void printPressure(ofstream &fid, TMesh* pTMesh_data){
	int nrows, ncols, i;
	pTMesh_data->allVertices.getSize(nrows,ncols);

	fid << "\nPOINT_DATA "<< nrows << endl;
	fid << "SCALARS Pressure float 1\n";
	fid << "LOOKUP_TABLE default\n";

	for(i=0;i<nrows;i++){
		VertexData* vdata = pTMesh_data->allVertices.getValue(i,0);
		fid << vdata->p << endl;
	}
}

void printSaturation(ofstream &fid, TMesh* pTMesh_data){
	// print saturation
	// = = = = = = = = = = = = = = = =  = = = = = = = = = = = = = = = = = =
	fid << "SCALARS Saturation float 1\n";
	fid << "LOOKUP_TABLE default\n";
	int nrows, ncols, i;
	pTMesh_data->allVertices.getSize(nrows,ncols);
	for(i=0;i<nrows;i++){
		VertexData* vdata = pTMesh_data->allVertices.getValue(i,0);
		fid << vdata->Sw << endl;
	}
}

//void printVelocity(ofstream &fid, TMesh* pTMesh_data){
//	fid << "VECTORS Velocity float 1\n";
//	int nrows, ncols, i, j;
//	int dim = pTMesh_data->getMeshDim();
//	double vel;
//	pTMesh_data->allVertices.getSize(nrows,ncols);
//	for(i=0;i<nrows;i++){
//		VertexData* vdata = pTMesh_data->allVertices.getValue(i,0);
//		for(j=0;j<dim;j++){
//			vel = vdata->vel[i];
//			if ( fabs(vel)<1e-8 ){
//				vel = .0;
//			}
//			fid << vel << " ";
//		}
//		fid << endl;
//	}
//}

void printPressureGradient(ofstream &fid, TMesh* pTMesh_data){
	fid << "VECTORS p_grad float\n";
	int nrows, ncols, i, j;
	pTMesh_data->allVertices.getSize(nrows,ncols);
	for(i=0;i<nrows;i++){
		VertexData* vdata = pTMesh_data->allVertices.getValue(i,0);
		double grad[3] = {vdata->grad_p[0],vdata->grad_p[1],.0};
		for(j=0;j<3;j++){
			fid << grad[j] << " ";
		}
		fid << endl;
	}
}

void printSaturationGradient(ofstream &fid, TMesh* pTMesh_data){
	fid << "VECTORS Sw_grad float\n";
	int nrows, ncols, i, j;
	pTMesh_data->allVertices.getSize(nrows,ncols);
	for(i=0;i<nrows;i++){
		VertexData* vdata = pTMesh_data->allVertices.getValue(i,0);
		double grad[3] = {vdata->grad_Sw[0],vdata->grad_Sw[1],.0};
		for(j=0;j<3;j++){
			fid << grad[j] << " ";
		}
		fid << endl;
	}
}

void printVelocity(ofstream &fid, TMesh* pTMesh_data){
	fid << "VECTORS Velocity float\n";
	int nrows, ncols, i, j;
	pTMesh_data->allVertices.getSize(nrows,ncols);
	for(i=0;i<nrows;i++){
		VertexData* vdata = pTMesh_data->allVertices.getValue(i,0);
		double grad[3] = {vdata->vel[0],vdata->vel[1],.0};
		for(j=0;j<3;j++){
			fid << grad[j] << " ";
		}
		fid << endl;
	}
}
